The Josephson effect throughout the BCS-BEC crossover 
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We study the stationary Josephson effect for neutral fermions across the BCS-BEC crossover, 
by solving numerically the Bogoliubov-de Gennes equations at zero temperature. The Josephson 
current is found to be considerably enhanced for all barriers at about unitarity. For vanishing 
barrier, the Josephson critical current approaches the Landau limiting value which, depending on 
the coupling, is determined by either pair-breaking or sound-mode excitations. In the coupling range 
from the BCS limit to unitarity, a procedure is proposed to extract the pairing gap from the Landau 
limiting current. 
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There has been much interest lately in the BCS-BEC 
crossover, motivated by its experimental realization with 
ultracold trapped Fermi atoms. Froperties that have 
been measured include radio-frequency spectra aimed at 
extracting, the pairing gap [1, 2, 3] and vortices under 
rotation [J]. In superconductors, these properties corre- 
spond to the occurrence of the gap parameter and the 
partial Meissner effect, which reveal the superconducting 
state. 

The two limiting situations, of fermions described by 
BCS theory in weak coupling and of composite bosons 
undergoing Bose-Einstein condensation in strong cou- 
pling, are connected with continuity because they cor- 
respond to the same spontaneous breaking of the gauge 
symmetry associated with the phase of the complex order 
parameter. In superconductors, this symmetry breaking 
is embodied by the Josephson effect, whereby a finite 
phase difference between two superconductors coupled 
via a link results in a steady current flow [^. 

Bosonic Josephson junctions have recently been real- 
ized with Bose-Einstein condensates Q. One can there- 
fore foresee an experimental study of the Josephson effect 
also with ultracold trapped Fermi atoms, whose mutual 
attraction can be varied from weak to strong coupling 
with the use of a Fano-Feshbach resonance [7|. In this 
context, the recent experimental study of stationary vs 
non-stationary flow of a bosonic condensate through a 
potential barrier may serve as a guide for a corre- 
sponding experiment with ultracold Fermi atoms, aimed 
at measuring specifically the Josephson effect. 

No theoretical study of the Josephson effect through- 
out the BCS-BEC crossover exists thus far, following the 
evolution between its Fermi and Bose versions. Pur- 
pose of the present work is to fill at least partially this 
gap, by calculating the Josephson current /phase relation 
throughout the BCS-BEC crossover at zero temperature 
within a mean-field approach for an inhomogeneous sys- 
tem in the presence of a barrier. The nontrivial results 
obtained by this study are, however, expected to hold 
even beyond mean field. These results regard specifi- 
cally the behavior of the Josephson characteristics and 



the maximum Josephson current for varying couplings 
and barriers, as well the relation of this maximum cur- 
rent to the Landau criterion for superfluidity across the 
BCS-BEC crossover. Our study will, in fact, reveal an in- 
timate connection between the maximum Josephson and 
Landau critical currents throughout the crossover. 

We consider a system of neutral fermions, mutu- 
ally interacting via an attractive short-range potential 
parametrized in terms of the scattering length ap- The 
BCS-BEC crossover is thus driven by the dimensionless 
coupling parameter {kFCtF)~^ where kp is the Fermi wave 
vector. This parameter ranges from (kpap)^^ <^ —1 in 
the weak-coupling (BCS) limit to 1 <C (kpap)^^ in the 
strong-coupling (BEC) limit, while the crossover region 
of interest is restricted to —1 < {kpap)~^ < -1-1 with the 
unitary limit at {kpap)~^ = 0. 

To be specific, we realize the Josephson link with a slab 
geometry, whereby a potential barrier V{x) of width L 
and height Vb > is embedded in a homogeneous super- 
fluid extending to inflnity on both sides of the barrier. Al- 
though the profile of the barrier is one-dimensional (with 
the X direction orthogonal to the barrier), the slab is fully 
three-dimensional as it extends along the two other (y 
and z) directions parallel to the barrier. This is because, 
on the BEC side of the crossover, the formation of com- 
posite bosons out of their fermionic constituents requires 
one to include wave vectors with components along all 
three dimensions. This marks a difference from previous 
treatments of the Josephson effect, that considered either 
fermions in the weak-coupling limit only M or point-like 
bosons whose formation is not an issue |lO'| . In addition, 
we regard the fermionic attraction to extend unmodified 
in the barrier region, thus anticipating the situation with 
ultracold Fermi atoms. 

To describe the effects of the potential barrier on the 
order parameter A(r) (where r = (x,y,z)), we solve 
the Bogoliubov-de Gennes (B-dG) equations for the two- 
component fermionic wave functions [lll |: 



H(r) 
A(r)' 



A(r) 
-H{r) 



Ui/(r) 



Ui/(r) 



(1) 



2 



Here, W(r) = — V^/(2m) + V{x) — fj, where m and fj. are 
the fermion mass and chemical potential. The function 
A(r) is determined via the self-consistency condition: 



A(r) 



(2) 



where —g is the strength of the local fermionic attraction, 
eliminated eventually in favor of the scattering length 
12l |. These equations have previously been consid- 



ered in the weak-coupling (BCS) limit, to determine the 
Josephson current for one-dimensional situations . 

As for any approach to a crossover problem, the present 
study of the Josephson effect has also to rely on a definite 
benchmark in the strong-coupling (BEC) limit. In this 
context, it has been shown that the fermionic B-dG 
equations ([1]) and ([2]) can be suitably mapped onto the 
Gross-Pitaevskii (GP) equation 14l llq 



-— + 2V{x) 
4m 
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$(r) = ^lB^{v) (3) 



for the condensate wave function $(r) of composite 
bosons with mass 2m, which form in the strong-coupling 
limit of the fermionic attraction. In the above equation, 
V{x) is the same potential of Eqs.lJ) and /x_b = 2/z-|-eo is 
the chemical potential for composite bosons (eq > be- 
ing the two- fermion binding energy such that fiB *C eo)- 
The two functions $(r) of Eq.([3l) and A(r) of Eq.® 
are related by <I>(r) = A(r)A/(m^aF)/87r. This mapping 
has been established for low-enough temperatures, so 
that all composite bosons reside in the condensate. This 
restriction matches one's expectation that a mean-field 
approach suffices to describe (at least qualitatively) the 
BCS-BEC crossover at low enough temperatures. 

The self-consistent solution of the B-dG equations first 
proceeds by extending to arbitrary values of the fermionic 
attraction the numerical approach introduced in Ref . [l6| 
for weak coupling. The profile of A(a;) is made piece- 
wise constant over a dense number of intervals (typi- 
cally 80), in such a way that the scattering problem is 
solved by elementary methods within each interval for 
given fermionic quasi-particle energy, while the scattering 
wave functions in contiguous intervals are connected by 
continuity conditions. Outgoing boundary conditions for 
waves impinging on the barrier from the left and right are 
used to identify the contribution of the continuos energies 
of Eqs.([T|). The contribution of bound levels (known as 
Andreev-Saint James states [17|) needs also be included. 
Self-consistency for the gap profile is implemented from 
the outcomes of such a scattering problem over a less 
dense grid of points (typically 20). Appropriate numer- 
ical strategies are required to speed up self-consistency, 
which is then reached with a limited number of cycles 
(typically 5). 

Numerical calculations have been performed by impos- 
ing either the value noq/m of the Josephson current (pro- 
vided it can be self-consistently sustained) where uq is 
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FIG. 1: (a) Magnitude 1A(3::)| and (b) phase 4i{x) of the order 
parameter vs the coordinate x orthogonal to the barrier with 
Lkp ~ 4.0 and Vo/Ef ~ 0.1, for the couplings {kpCLF)"^: 
—0.8 (dashed line); 0.0 (dotted Une); and 3.0 (dash-dotted 
line). The last curve is compared with the independent solu- 
tion of the GP equation (full line). Here Ep = kp/{2m) is 
the Fermi energy. The inset shows the density profiles n{x). 



the bulk density, or the value of the asymptotic phase 
difference (5</> = (f>{x = +oo) — (j){x = —oo) between the 
two sides of the barrier. Under these circumstances, we 
write A{x) = |A(a;)| exp[2iqx + i(j){x)]. 

Figures 1(a) and 1(b) present the spatial variation of 
the magnitude | A(a::)| and phase (/'(x) of the order param- 
eter across a barrier centered aX x — 0, for three char- 
acteristic values of (kpap)^^ and for the corresponding 
maximum value of the Josephson current. The indepen- 
dent solution of the GP equation is also compared with 
the results of the B-dG equations for (kpap)^^ = 3, when 
one expects the BEC limit to have been reached. One 
sees that the depression of |A(a:)| from its bulk value 
|Ao| due to the presence of the repulsive barrier causes 
a sharp variation of the phase (t>{x) across the barrier in 
order to keep the current uniform, thus resulting in the 
total phase difference 6(j) accumulating across the barrier. 
Both the depression of |A(a:)| at a; = and the value of 
6(j) systematically increase from weak to strong coupling. 
These features are reflected in the corresponding density 
profiles n{x) shown in the inset of Fig. 1(b). In the BCS 
limit, Friedel oscillations modulated by 2kp affect n{x) 
as well as |A(x)| and <j>{x). 

The characteristic Josephson current/phase relation 
J{5<j)) for a range of different coupling encompassing uni- 
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FIG. 2: (a) J{5(j>) for the same values of L and Vo of Fig.l 
and for the couphngs (/cfaf)^^: —0.8 (squares); 0.0 (circles); 
1.0 (triangles), (b) J{S(j}) at unitarity for Lkp = 4.0 and 
different values of Vb/i?F: 0.025 (squares); 0.10 (circles); 0.40 
(triangles). [Dashed lines interpolate among data points.] 



tarity is shown in Fig. 2 (a), for the same values of L and 
Vo considered in Fig.l. Note that the Josephson current 
is considerably enhanced at unitarity with respect to the 
BCS and BEC sides, and that the curves J{5(j)) stretch 
from being proportional to sin{S<p) when (fci?aF)~^ = 
1.0 to being (almost) proportional to cos((50/2) when 
(kpap)^^ — —0.8. By scanning the coupling parame- 
ter for given barrier, we have further verified that: (i) In 
the BEC limit one invariably obtains J{S(j)) oc sin((50), re- 
flecting the fact that in this limit the barrier heigth Vq be- 
comes larger than the relevant energy scale fj,B', (ii) In the 
BCS limit (whereby the coherence lenght by far exceeds 
the barrier width) one reduces to considering a delta- 
like barrier with strength Z ~ {Lkp){VQ/Ep), for which 
J{S(j)) oc sin((50) when Z ^ 1 and J{S(j)) oc cos((50/2) 
when Z <C 1. Such a progressive evolution from a 
sin((50)-dependence for strong barriers to a cos((5(/)/2)- 
dependence for weak barriers occurs actually for any cou- 
pling, as shown for instance at unitarity in Fig. 2(b) [l^ . 

From these Josephson characteristics one can extract 
the maximum (critical) current Jc = noqc/m for given 
barrier and coupling. The corresponding critical veloc- 
ity qc/m is shown in Fig. 3 vs {kpaF)~^ from weak to 
strong coupling and for several barrier heights. All curves 
present a maximum at about unitarity, with the value of 
the maximum increasing as the ratio Vo/Ep is decreased. 
The question thus naturally arises whether this limiting 
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FIG. 3: Maximum velocity qc/m vs {kpap)'^ for different 
barriers with Lkp — 4.0 and Vo/Ep: 0.025 (squares); 0.10 
(circles); 0.20 (stars); 0.40 (triangles). [Dashed lines interpo- 
late among data points.] Full lines represent the two Landau 
branches for pair breaking (left) and sound mode (right). 



procedure of lowering the barrier height leads to an in- 
trinsic upper value of qc/m for given coupling. 

On physical grounds, a vanishingly small barrier acts 
as an impurity that probes the stability of the homoge- 
neous superfluid flow, thus playing a similar role to the 
walls of the container in the context of the Landau crite- 
rion for superfluidity [l^ . One then expects the velocities 
reported in Fig. 3 to never exceed the corresponding value 
of the critical velocity obtained by the Landau criterion. 
This critical velocity is determined by the available quasi- 
particle excitations which are of a different nature on the 
two sides of the crossover. Specifically, on the BCS side 
there occur the pair-breaking excitations characteristic of 
BCS theory, while on the BEC side one expects sound- 
mode quanta to be the relevant excitations. 

Accordingly, in Fig. 3 we have also plotted (full lines) 
the critical velocities corresponding to the above two 
branches of the Landau criterion. In particular, the 
"left" branch corresponds to the expression q^/w = 
y/ fi'^ + Aq — fj,, which is obtained by finding the small- 
est value of qc for which E{k) — {qc/m)k — has solu- 
tions, where E(k) = ^(fc2/(2TO) - ^)2 is the BCS 
single-particle energy of wave vector k when the super- 
fluid is at rest. The "right" branch corresponds instead 
to the sound velocity of the Bogohubov- Anderson mode 
obtained within the BCS-RPA approximation [i^l , which 
generalizes the result q^/Tn = hb valid in the BEC limit. 
These two Landau branches represent an upper limit for 
the qc/m vs lkpap)~^ curves obtained from our numer- 
ical calculations in the presence of a barrier, the latter 
consistently tending to the former from below as the bar- 
rier height is progressively decreased. The two Landau 
branches intersect each other at about unitarity, where 
the maximum value of qc/m is then achieved 21, 2^ . 

A few comments on how these two Landau branches 
result from the solution of the B-dG equations are in 
order. As a consequence of Galilean invariance, no up- 
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per critical velocity would result from the equations ^ 
and ([2]) for a strictly homogeneous system. It is thus only 
through a limiting procedure, implemented, e.g., by a pro- 
gressively vanishing barrier height, that an upper value of 
the critical velocity results when solving these equations. 
Yet, by this procedure one would naively expect to re- 
cover only the left Landau branch corresponding to pair- 
breaking excitations, arguing that the B-dG equations 
correspond to the BCS mean field which contains explic- 
itly only this type of excitations. The fact that the right 
Landau branch for sound-mode excitations emerges, too, 
from our numerical calculations appears thus remarkable. 
The crucial point is that, in the presence of a non-trivial 
geometry, the imprint of the excitation spectrum is found 
already in the ground-state wave function at given (su- 
per) current, which in the present context corresponds 
to the mean-field equations ([T]) and ([2]) at zero temper- 
ature. Since these equations reduce to the GP equation 
in the EEC limit, also their excitations on the BEG side 
correctly correspond to those of the GP equation [i^ . 

The above results suggest us a procedure to extract the 
value of the pairing gap (order parameter) from the BCS 
limit to unitarity, which is the re gion of most interest 
for this quantity. Previous studies [2^ have shown that 
for pair-breaking excitations the BGS expression E{k) 
remains valid even beyond mean-field, with appropriate 
values of Aq and /i. In this way, the Landau expres- 
sion for Qc reported above remains also valid, and can be 
inverted to give Ao/Ep — 2[qc/kF)^ {qc/kpY + [ijEp. 
The value of At^jEp is thus obtained from the exper- 
imental value of qc/kp (resulting, for instance, by ex- 
tending to ultracold Fermi atoms the method recently 
used to study the flow of a bosonic condensate through a 
potential barrier [s*]) and the value of ^i/Ep determined 
either experimentally or from theoretical calculations. In 
this respect, we have verified that the differences in the 
values of Aq obtained by using mean-field or QMC val- 
ues of /i never exceed 10% over the relevant range of 
the coupling parameter. In addition, the Landau value 
of qc is expected to be measured with reasonable accu- 
racy on the BGS side of the crossover, where this limiting 
value is reached rather quickly by decreasing Vq as shown 
in Fig. 3. Measurement of qc would therefore provide a 
quantitative probe for fermionic superfiuidity. 

We are indebted to G. Deutscher, R. Hulet, and F. 
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[2] M. Greiner, C.A. Regal, and D.S. Jin, Phys. Rev. Lett. 

94, 070403 (2005). 
[3] C.H. Schunck, Y Shin, A . Schirotzek, M.W. Zwierlein, 

and W. Ketterle, 'cond-mat/0702066 
[4] M.W. Zwierlein, J.R. Abo-Shaeer, A. Schirotzek, C.H. 

Schunck, and W. Ketterle, Nature 435, 1047 (2005). 
[5] B.D. Josephson, in Superconductivity, R.D. Parks, Ed. 

(Dekker, New York, 1969), p.423. 
[6] M. Albiez et al, Phys. Rev. Lett. 95, 010402 (2005); see 

also F.S. Cataliotti et al. Science 293, 843 (2001). 
[7] S. Inouye et al. Nature 392, 151 (1998). 
[8] P. Engels and C. Atherton, arXiv;0704.2427vl[cond- 

mat. other]. 

[9] A. Furusaki and M. Tsukada, Phys. Rev. B 43, 10164 

(1991). See also G. Wendin and V.S. Shumeiko, Phys. 

Rev. B 53, R6006 (1996), and references quoted therein. 
[10] S. Raghavan, A. Smerzi, S. Fantoni, and S.R. Shenoy, 

Phys. Rev. A 59, 620 (1999). 
[11] P. G. De Gennes, Superconductivity of Metals and Alloys 

(Benjamin, New York, 1966). 
[12] C.A.R. Sa de Melo, M. Randeria, and J. R. Engelbrecht, 

Phys. Rev. Lett. 71, 3202 (1993). 
[13] P. Fieri and G. C. Strinati, Phys. Rev. Lett. 91, 030401 

(2003) . 

[14] E. P. Gross, II Nuovo Cimento 20, 454 (1961). 

[15] L. P. Pitaevskii, Sov. Phys. JETP 13, 451 (1961). 

[16] R.A. Riedel, Li-Fu Chang, and P.P. BagweU, Phys. Rev. 
B 54, 16082 (1996). 

[17] A. Andreev, Soviet Phys. JETP 19, 1228 (1964); D. 
Saint-James, J. de Physique 25, 899 (1964). 

[18] This evolution has been previously noted in the weak- 
coupling limit by F. Sols and J. Ferrer, Phys. Rev. B 49, 
15913 (1994). 

[19] See, e.g., A. A. Abrikosov, L.P. Gorkov, and I.E. 
Dzyaloshinski, Methods of Quantum Field Theory in Sta- 
tistical Physics (Dover, New York, 1963). 

[20] See M. Marini, F. Pistolesi, and G.C. Strinati, Eur. Phys. 
J. B 1, 151 (1998) for analytic expressions. 

[21] That superfluidity is most robust against suppressing 
mechanisms in the crossover region was pointed out by 
N. Andrenacci, P. Fieri, and G.C. Strinati, Phys. Rev. B 
68, 144507 (2003). See also R. Sensarma, M. Randeria, 
and T.-L. Ho, Phys. Rev. Lett. 96, 090403 (2006). 

[22] Taking into account the appropriate reduction of the 
bosonic scattering length from 2aF to 0.6aF [D.S. Petrov, 
C. Salomon, and G.V. Shlyapnikov, Phys. Rev. Lett. 93, 
090404 (2004)] will only quantitatively modify the right- 
hand side of the curves of Fig. 3, while preserving the 
physical conclusions of the present calculation. 

[23] G.C. Strinati and P. Fieri, Phys. Rev. A 69, 011601(R) 

(2004) . 

[24] P. Fieri, L. Pisani, and G.C. Strinati, Phys. Rev. B 70, 
094508 (2004); S.Y. Chang, V.R. Pandharipande, J. Carl- 
son, and K.E. Schmidt, Phys. Rev. A 70, 043602 (2004). 



[1] C. Chin et al, Science 305, 1128 (2004). 



